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1.  Technical 

The  initial  focus  of  the  phase  two  effort  is  towards  developing  methods,  tools  and  databases  to 
support  the  various  environmental  data  required  by  the  RF  propagation  prediction  models. 
Another  aspect  is  to  build  upon  the  Phase  I  accomplishments  to  develop  interfaces  that  allow 
different  models  to  be  coupled  together  in  order  to  do  mixed  path  propagation  computations. 
Lastly,  in  order  for  a  mixed  path  propagation  model  to  be  truly  useful  for  many  of  today’s 
problems,  methods  must  be  developed  to  couple  the  energy  radiated  from  antennas  to 
propagation  models.  This  is  especially  necessary  as  the  link  budget  equation  becomes  no  longer 
applicable. 

1.1.  Technical  Objectives  Accomplished 

1.1.1.  Task  Progress 

Over  the  past  few  months  it  has  become  increasingly  clear  that  the  customers  for  the  Mixed  Path 
project  desire  API  to  the  various  Wireless  InSite  calculation  engines.  This  is  because  then  the  WI 
engines  can  be  incorporated  into  too  such  as  OPNET  and  QualNet  for  network  analysis  or  into 
various  signal  processing  software.  For  the  network  analysis  all  that  is  required  is  an  extremely 
fast  computing  propagation  model  that  captures  the  mixed  path  propagation  at  an  acceptable 
resolution  and  fidelity.  For  signal  processing  analyses,  the  propagation  channel  is  required.  Thus 
a  time  domain  capability  API  is  required  within  Wl.  Lastly  an  easy  to  use  graphics  GUI  is 
required  by  all  of  WI  customers  to  allow  the  scene  to  be  more  easily  rendered  and  modified. 
Further  Wireless  InSite  customers  require  the  ability  to  run  many  different  permutations  and 
combinations  of  scenarios.  Thus  scripting  is  being  developed  for  Wireless  InSite.  Thus  this 
report  shall  focus  on  accomplishments  in  these  areas. 

1.1.1. 1.  Fast  Urban  API 

1.1. 1.1.1.  Background 

In  order  for  Wireless  InSite  to  interface  to  network  simulations  such  as  OPNET  and  QualNet,  it 
is  necessary  for  an  API  to  be  constructed  for  the  propagation  models.  An  API  could  be  made 
directly  to  point  to  point  version  of  the  fast  propagation  models.  The  calculation  API  originally 
conceived  for  newly  developed  faster  ray  based  methods  to  support  rapid  computation. 

1.1. 1.1.2.  Status 

This  capability  is  to  be  available  as  a  commercial  product  in  August  2006. 

1.1. 1.1.3.  Technical  Details 

See  Section  4.2 

1.1. 1.2.  Faster  urban  propagation  models 
1.1. 1.2.1.  Background 

In  order  to  support  rapid  urban  calculations  suitable  for  integration  into  networks  simulations,  a 
rapid,  perhaps  lower  fidelity,  urban  propagation  models  are  required.  This  model  shall  be  site 
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specific.  Remcom’s  high  fidelity  models  found  in  Wireless  InSite  could  be  used  to  form  the  basis 
of  semi-empirical  models  or  be  used  for  validation. 

1.1. 1.2.2.  Current  Status 

Remcom’s  vertical  plane  urban  propagation  model  is  complete  and  integrated  into  the  fast  urban 
API. 


1.1. 1.2.3.  Technical  Details 

See  Section  4. 1 . 

1.1. 1.3.  Two  Dimensional  Building  Simplification 

1.1. 1.3.1.  Background 

While  Wireless  InSite  currently  has  a  two  dimensional  building  simplification  algorithm,  it 
requires  extensive  user  interaction  and  proficiency.  Further  the  current  simplification  metrics 
must  be  applied  one  building  at  a  time.  This  makes  building  simplification  arduous  and  time 
consuming.  This  task  is  to  find  better  methods  that  require  little  user  interaction  and  can  be 
applied  over  the  entire  area  of  interest. 

1.1. 1.3.2.  Current  Status 

This  was  integrated  into  Wireless  InSite. 

1.1. 1.3.3.  Technical  Details 

See  Section  4.3.1. 

1.1. 1.4.  Three  Dimensional  Building  Simplification 

1.1. 1.4.1.  Three  Dimensional  Simplification  Background 

Often  the  buildings  in  an  urban  or  suburban  scenario  are  too  small/short  for  the  urban  canyon’s 
propagation  methods  to  function  properly.  Further  time  consuming  3-dimensional  ray  tracing 
computations  are  directly  related  to  the  number  of  faces  found  in  the  scenario.  Thus  it  is  desired 
to  find  methods  to  reduce  the  number  of  faces  in  the  scenario  to  a  reasonable  yet  representative 
realization. 

1.1. 1.4.2.  Three  Dimensional  Simplification  Background  Status 

Because  of  the  technical  complexity,  this  task  was  discontinued. 

1.1. 1.4.3.  Technical  Details 

See  Section  4.3.2. 

1.1. 1.5.  Soil  Dielectrics 
1.1. 1.5.1.  Background 

In  order  to  support  an  easy  method  to  determine  the  world  wide  ground  dielectric  properties, 
Remcom  investigated  possible  methods  to  compute  the  dielectric  constants  from  the  soil 
constituents  and  moisture  content. 


2 


Prepared  by  Remcom,  Inc. 


Final  Report 


1.1. 1.5.2.  Current  Status 

The  algorithms  are  complete.  Software  to  determine  the  opil  constituents  worldwide  is  complete. 

1.1. 1.5.3.  Technical  Details 

See  Section  4.4. 

1.1. 1.6.  Wireless  InSite  Optimization 

The  goal  is  to  significantly  reduce  the  computation  time  required  by  Wireless  InSite’s  full  three- 
dimensional  (3D)  propagation  model.  This  will  be  accomplished  by  decreasing  the  run  time  of 
the  SBR  ray  tracing  procedure.  Other  goals  include  creating  a  calculation  engine  executable  for 
a  32-bit  Linux  operating  system,  and  to  enable  parallel  processing  on  a  Linux  cluster. 

1.1. 1.6.1.  Reduce  Run  Time  of  Full  3D  Propagation  Model 

Tables  1  and  2  in  Section  1.2.6  shows  the  speedup  for  version  2.2,  and  for  the  version  of  2.3  that 
we  are  currently  testing.  The  speedup  for  the  2.2  version  is  usually  in  the  5-9  range  and  generally 
improves  slightly  as  the  number  of  buildings  increases.  The  speedup  for  version  2.3  is  currently 
in  the  7-13  range.  We  hope  to  improve  that  a  little  more  before  we  release  version  2.3  at  the  end 
of  May. 

1.1. 1.6.2.  Shapeflie  Import 

Although  this  was  no  included  in  the  current  contract,  we  have  decided  to  add  the  ability  to 
version  2.3  to  import  building  data  from  shapefiles  which  use  the  polygon  and  polygon  Z  shape 
types.  We  will  import  the  data  directly  from  the  shapeflie,  not  from  the  intermediate  DXF  file 
that  you  needed  to  use  with  version  2.1.  We  will  support  files  with  building  coordinates  in 
longitude/latitude,  UTM  and  local  Cartesian  frames  that  are  not  geo-referenced.  We  are  still 
investigating  whether  we  will  be  able  to  import  data  using  other  map  projections.  Version  2.2 
does  not  have  any  capabilities  for  directly  importing  shapefiles. 

1.1. 1.7.  Calculation  Engine  for  Linux  Operating  Systems 

The  2.3  calculation  engine  will  run  on  a  computer  with  Linux  operating  system.  The  software 
can  be  licensed  using  either  a  node  locked  license  or  a  floating  license.  The  calculation  engine 
can  be  executed  either  from  the  command  line  or  by  means  of  a  script  file.  The  name  of  the 
project  file  to  run  is  one  of  the  command  line  options  that  the  user  can  specify.  The  GUI  will  still 
only  run  on  a  Windows  operating  system. 

1.1. 1.7.1.  Running  Calculations  on  Linux  Clusters 

The  2.3  Wireless  InSite  calculation  engine  now  runs  on  a  Linux  cluster.  The  GUI  will  generate 
script  files  for  submitting  jobs  to  a  cluster  using  PBS  scheduler.  All  output  can  be  viewed  using 
the  GUI  once  the  files  are  moved  to  a  location  which  the  GLH  can  access. 

1.1. 1.8.  Building  Pre-Processor 

We  have  made  a  number  of  improvements  to  the  building  preprocessor  (building  simplifier). 
Many  of  the  new  capabilities  are  available  in  version  2.2 
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1.1. 1.9.  Test  Cases  and  Run  Time  Benchmarks 

We  have  developed  several  test  cases  for  comparing  the  new  run  times  to  the  2.1  run  times  in 
urban  environments.  The  tests  cases  use  a  variety  of  building  shapes,  transmitter  heights  and 
locations,  and  numbers  of  receivers.  The  test  cases  are  also  used  to  verify  that  the  run  tim.e 
reductions  are  not  affecting  the  accuracy  of  the  results.  Run  times  for  downtown  Chicago  test 
cases  are  discussed  below. 

1.1. 1.9.1.  Description 

The  section  of  downtown  Chicago  used  for  the  run  time  testing  is  shown  in  Figures  1  and  2. 
This  area  is  about  1  km  x  1  km  and  contains  about  100  buildings  with  2200  facets.  The  locations 
of  the  five  transmitter  sites  are  indicated  by  the  green  boxes.  The  calculations  were  performed 
for  1850  receiver  locations  spaced  at  5  meters  along  the  main  streets  in  this  area  with  only  one 
transmitter  active  at  a  time.  The  path  loss  predictions  along  the  streets  for  each  transmitter 
location  are  shown  in  Figures  3-7. 

The  run  times  for  each  transmitter  location  were  recorded  using  versions  2.1  and  2.2,  and  the 
preliminary  version  of  2.3.  The  run  times  and  the  computed  speedups  for  versions  2.2  and  2.3 
relative  to  2.1  are  shown  in  Table  1  and  2.  The  average  speedup  for  version  2.2  relative  to 
version  2.1  is  about  7.0  for  the  total  run  time,  and  about  8.0  for  just  the  ray-tracing  portion  of  the 
calculation.  The  average  speedup  for  version  2.3  relative  to  version  2.1  is  about  10.4  for  the 
total  time,  and  about  1 1.4  for  just  the  ray  tracing. 
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Project  view;  (Do v/n town  Chicago) 


Figure  1:  Plan  view  of  downtown  Chicago  area  with  the  transmitter  and  receiver  locations 


Eroject  Edt  View  Output  Select 
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Figure  2:  Three-dimensional  view  of  downtown  Chicago  area  with  transmitter  and  receiver 

locations 


1.1. 1.9.2.  Runtime  and  Speedup 


1.1. 1.9.3.  Version  2.2  Run  Time  and  Speedup 

Table  1:  Run  times  and  speedup  for  version  2.2  on  Downtown  Chicago  test  cases 


Test  Case 

Tota 

Run  Time  ( 

li:m:s) 

Ray  Tracing  Run  Time  (h:m:s) 

v2.1 

v2.2 

Speedup 

v2.1* 

v2.2 

Speedup* 

Street  Level  Tx  #1 

8:57:42 

1:28:42 

6.1 

8:52:05 

1:24:03 

6.3 

Street  Level  Tx  #2 

14:14:00 

2:25:06 

5.9 

14:06:42 

2:13:15 

6.4 

Rooftop  Tx  #1 

10:32:00 

1:26:52 

7.3 

10:22:38 

1:14:44 

8.3 

Rooftop  Tx  #2 

8:31:00 

1:05:19 

7.8 

8:24:01 

0:56:15 

9.0 

High  Tower  Tx  #1 

13:58:00 

1:39:56 

8.4 

13:54:07 

1:36:47 

8.6 

*  Estimate 
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1.1. 1.9.4.  Version  2.3a  Run  Time  and  Speedup 


Table  2:  Run  times  and  speedup  for  version  2.3a  on  Downtown  Chicago  test  cases 


Test  Case 

Tota 

Run  Time  l 

h:m:s) 

Ray  Tracing  Run  Time  (h:m:s)  | 

v2.1 

v2.3a 

Speedup 

v2.1* 

v2.3a 

Street  Level  Tx  #1 

8:57:42 

1:04:12 

8.4 

8:52:05 

1:01:01 

8.7 

Street  Level  Tx  #2 

14:14:00 

1:34:44 

9.0 

14:06:42 

1:23:50 

10.1 

10:32:00 

0:59:00 

10.7 

10:22:38 

0:47:37 

11.4 

8:31:00 

0:45:16 

11.3 

8:24:01 

0:36:39 

13.1 

13:58:00 

1:05:48 

12.7 

13:54:07 

1:00:00 

13.9 

*  Estimate 
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1.1. 1.9.5.  Path  Loss  Predictions 


Figure  3:  Path  loss  predictions  along  streets  for  Street  Level  Transmitter  #1  at  600  MHz 
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Project  view:  (Downtown  Cfikago) 


Ready. 


Figure  4:  Path  loss  predictions  along  streets  for  Street  Level  Transmitter  #2  at  600  MHz 
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Project  view:  (Downtown  Chicago) 


Eroject  Efit  yiew  QMtput  Seiect  bMp 


-160  0ijB 


Figure  5:  Path  loss  predictions  along  streets  for  Rooftop  Transmitter  #1  at  600  MHz 
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If^  Project  view:  (Downtown  CNcago) 


Eroject  E^lt  y)ew  Qiutput  Select  tldp 
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Figure  6:  Path  loss  predictions  along  streets  for  Rooftop  Transmitter  #2  at  600  MHz 
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Figure  7:  Path  loss  predictions  along  streets  for  the  High  Tower  Transmitter  at  600  MHz 


The  goal  is  to  significantly  reduce  the  computation  time  required  by  Wireless  InSite’s  full  three- 
dimensional  (3D)  propagation  model.  This  will  be  accomplished  by  decreasing  the  run  time  of 
the  SBR  ray  tracing  procedure.  Other  goals  include  creating  a  calculation  engine  executable  for 
a  32-bit  Linux  operating  system,  and  to  enable  parallel  processing  on  a  Linux  cluster. 

1.1. 1.9.6.  Timing  Functions 

During  this  month  we  began  work  on  creating  “timer”  functions  that  can  be  inserted  around  any 
section  of  code.  These  functions  are  inserted  permanently  into  the  source  code  for  accurately 
monitoring  run  time  of  different  components.  A  summary  of  all  run  time  data  this  is  collected 
will  be  written  to  the  diagnostics  (.diag)  file  when  the  calculation  finishes.  When  the  calculation 
engine  is  run  from  the  GUI,  the  run  time  diagnostics  will  also  be  written  to  the  window  which 
displays  calculation  progress  (calculation  log). 
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The  C-H-  language  was  used  to  develop  classes  capable  of  timing  sections  of  code  and  of 
reporting  the  collected  timings  in  a  unified  format.  An  interface  was  developed  between  the  C++ 
classes  and  the  Fortran  code  used  in  the  calculation  engine.  103  important  locations  in  the 
calculation  engine  were  identified  and  a  timer  was  placed  in  each  location.  InSite  projects  were 
built  that  use  parts  of  the  calculation  engine  that  AFIWC  will  use.  Those  projects  were  run  and 
the  timing  results  were  analyzed. 

The  analysis  successfully  identified  the  critical  path  of  the  computation  and  sections  of  the 
calculation  engine  important  to  the  overall  run  time.  Critical  sections  were  closely  examined. 
Small  changes  in  those  sections,  most  of  which  were  in  the  ray-tracing  portion  of  the  engine, 
proved  effective  in  reducing  total  run  time.  The  extent  of  improvement  varies  by  project.  Some 
projects  currently  run  17-22%  faster.  A  developer  reports  that  one  project  runs  in  25%  of  the 
original  time,  now  running  in  15  minutes  compared  to  previously  running  in  an  hour.  Projects 
with  more  diffractions  appear  to  have  more  improvement.  Note  that  some  portion  of  the 
improvements  come  from  using  a  new  version  of  the  Lahey  Fortran  compiler. 

Summary  timing  data  now  is  written  to  the  calculation  log  window,  visible  in  the  graphical  user 
interface  whenever  a  calculation  is  running.  Timings  are  shown  with  each  percentage  complete 
indicator  written  to  the  calculation  log  window.  Detailed  timing  data  is  also  written  in  the  .diag 
file  associated  with  each  project. 

The  work  on  this  task  has  been  completed. 

1.1.1.10.  Develop  a  Faster  3D  SBR  Ray-Tracing  Algorithm 

Because  of  the  similarities  of  the  SBR  approach  used  in  WI  to  the  ray  tracing  performed  for 
computer  graphics  applications  involving  the  rendering  of  complex  scenes,  there  are  a  number  of 
promising  techniques  available  for  deceasing  the  computation  time.  We  have  begun  to  adapt 
recent  developments  in  computer  graphics  techniques  to  the  ray  tracing  required  by 
electromagnetics  computations.  The  work  on  this  task  will  continue  into  next  month. 

With  the  integration  of  C++  in  the  calculation  engine  described  in  section  1.1. 1.9. 6,  the 
foundation  has  been  laid  for  replacing  the  ray-tracing  portion  of  the  Fortran  engine  with 
Remcom's  highly-optimized  C++  Optimized  Ray-Tracing  Engine  (ORE).  A  design  for 
integrating  ORE  has  been  reviewed  and  approved  by  the  InSite  Project  Leader. 

The  work  on  this  task  will  continue  into  next  month. 

1.1.1.10.1.  Calculation  Engine  Executable  for  Linux  Operating  Systems 

The  Linux  port  of  the  calculation  engine  is  in  progress.  A  set  of  regression  tests  has  been  put  in 
place  and  results  from  44  separate  runs  of  the  Linux  version  of  the  engine  are  regularly  compared 
to  results  from  the  Windows  version.  The  Linux  test  cases  are  run  on  Remcom’s  Linux  cluster 
under  PBS,  with  results  brought  back  to  the  Windows  computer  for  comparison  in  the  regression 
tests. 

The  test  cases  include  both  optimized  runs  and  debug  runs.  The  debug  runs  use  the  full  checking 
capabilities  of  the  Lahey  Fortran  compiler  to  ensure,  for  example,  that  array  indices  remain  in 
bounds  and  that  inputs  to  intrinsic  numeric  functions  are  within  the  correct  range. 
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The  work  on  this  task  will  continue  into  next  month. 

1.1.1.10.2.  Enhancements  to  Building  Pre-Processor 

The  building  pre-processor  is  used  to  eliminate  unnecessary  detail  from  urban  building  data. 
During  this  month  we  have  begun  testing  the  new  building  pre-processing  algorithms  and 
interface.  The  work  on  this  task  will  continue  into  next  month. 

The  work  on  this  task  will  continue  into  next  month. 

1.1.1.10.3.  Test  Cases 

During  this  month  we  have  modified  the  test  cases  which  will  be  used  in  the  future  to 
demonstrate  run  time  reduction.  Test  cases  have  been  developed  for  the  following  scenarios: 

a.  Dense  urban  area  with  street  level  transmitters  and  street  level  receivers 

b.  Dense  urban  area  with  rooftop  transmitters  and  street  level  receivers 

c.  Dense  urban  area  with  airborne  transmitters  above  rooftops  and  street  level  receivers. 

The  work  on  this  task  will  continue  into  next  month. 

1.2.  Immediate  Plans 

Begin  work  on  Wireless  2.3  capabilities. 

1.3.  Technical  Prc^blems 

N/A 

1.3.1.  Technical  Problems  Encountered 

N/A 

1.3.2.  Cost/Schedule  Impact 

N/A 

1.3.3.  Steps  taken  to  Remedy  Problems 

N/A 

1.4.  Resolutions  to  Previous  Technical  Problems 

N/A 

1.5.  Recommendations 

N/A 
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2.  Cost  Incurred  and  Outstanding  Commitments 

2.1.  Current  Costs  Incurred 

$0.00 

2.2.  Cumulative  Costs  Incurred 

$730,759.39 

2.3.  Funds  Remaining 

$15,801 

2.4.  Spending  Costs 

The  total  available  funding  on  this  project  is  $746,560.  Based  on  this  expectation,  we  have 
revised  our  spending  plan  as  shown  in  the  following  chart. 


Mixed  Path 


Time 
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3.  Project  Management  Review  (PMR) 

3.1.  Agendas 

N/A 

3.2.  Presentation  Materials 

N/A 

3.3.  Minutes 

N/A 

4.  Appendices  Technical  Reports 

4.1.  Vertical  Plane  Urban  Propagation  Model 

The  essential  characteristics  of  the  empirical  models  originally  attributable  to  Bertoni  and 
Ikegami  are  identified  and  incorporated  as  components  of  a  more  general  model  which  provides 
alternate  formulations  of  the  established  algorithms  for  multi-screen  diffraction  (MSD),  roof-to- 
street  attenuation  and  novel  formulations  for  base  station  elevations  below  that  of  roof-top  level. 
The  applicability  of  the  conventional  elements  is  enhanced  with  supplementary  methods 
appropriate  for  more  general  urban  terrain  not  confined  to  collinear  buildings  of  uniform  spacing, 
thus  permitting  more  general  ray  paths  composed  of  isolated  MSD  diffraction  configurations  and 
conventional  (optically  independent)  multiple  knife  edge  diffracting  elements. 

4.1.1.  Modifications  and  Extensions  of  Conventional  Methods  of  Vertical 
Plane  Urban  Propagation 

4.1. 1.1.  Alternative  Formulation  of  Multi-Screen  Diffraction  Loss  for  Base 
Stations  of  Arbitrary  Elevation  Relative  to  Neighboring  Buildings 

The  Walisch  Ikegami  model  [6]  is  valid  for  urban  terrain  conforming  to  the  assumptions  of 
uniform  spacing  and  of  nearly  collinear  roof  tops  as  in  Figure  8.  The  model  is  unreliable  for 
transmitter  elevation  approaching  roof  top  level  and  performs  poorly  for  base  station  heights 
below  roof  top  level.  Acknowledging  these  deficiencies,  the  proposed  model  provides  an 
accurate  and  numerically  tractable  alternative  formulation  of  multi-screen  diffraction  (MSD) 
based  on  the  exact  solutions  of  Boesma  [1]  and  valid  at  arbitrary  transmitter  elevations  relative  to 
the  line  of  screens  used  to  approximate  buildings  in  an  urban  terrain.  The  proposed  alternate 
formulation  of  the  MSD  loss  is  expressed  as  a  height  dependent  correction  of  the  value  for 
collinear  transmitter  and  roof  tops  and  is  of  the  explicit  form. 
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Figure  8 


a  "'(gj-  Qc  ""(0) + /<  [a  )  -  0c  ""(0)],  (g.  <  0) 


(1) 


where 


0c  ""(^c  )  -  0c  (0)  +  />  10c  (go  )  -  0c  *"(0)],  (gc  >  0) 

l0g(l+ I  gc  I) 


/>  = 


/<  = 


log(l+  I  go  I) 


log(l+|gc  I) 


(2) 

(3) 

(4) 


log(l+ 1  I) 

and  where  =2  and  where  g^  is  the  (scaled)  deviation  (Ay)  of  the  transmitter  height 

from  roof  top  level  ,  g^  =Ay/yl{Ad)  ,  Q^'^iO)  is  the  exact  loss  (-20 log, ^ (1/M))  for  incidence 
upon  the  edge  of  the  Mth  screen  for  a  transmitter  aligned  exactly  with  the  (M)  roofs  and 
0c*^(gc<)(0c*^(go))  ^  numerical  approximation  of  the  exact  value  obtained  using  Boesma’s 

formulation  as  evaluated  for  g^  -gc<(gc>)-  addition,  the  parameter  p  appearing  in  Eq.(2)  is 
computed  according  to  the  relation 

P.(gc)  =  ^^;(gc<(i/2)) 
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P.(^.)  =  0.8;(g,>(l/2)) 


The  two  bounding  interpolation  control  quantities  Qc'^igo)^  3nd  which  are  the 

values  of  loss  at  two  reference  values  of  the  scaled  arguments  and  ,  are  independently 
evaluated  according  to  a  similar  approximate  expressions. 

a  "(«„ ) = a  "(g.« ) + la  )  -  a "(«,« )i  (s) 

a"te„)=a"(&>.)+(i>>g([^!jj^^^]))'[a"(g„.)-a"te,..)i  «;) 

where  q  is  an  empirical  parameter  optimized  as  q  =  \.l .  Other  optimized  numerical  parameters 
required  in  the  approximation  of  Boesma’s  expressions  include 

e/(^e«)  =  19.0909  =  60.5257, =  -0.74775,  =  22.52014  (7) 

Specific  results  for  the  approximate  expression  and  the  exact  equivalent  due  to  Boesma  appear  in 
Figure  9  in  which  the  multi-screen  diffraction  loss  is  plotted  as  a  function  of  the  scaled 

height  displacement  coordinate  for  values  of  the  number  of  screens  (knife  edges)  in  multiples 
of  two  between  two  and  100.  It  is  noted  that  for  the  minimum  number  of  screens  considered, 
M=2,  the  physical  interpretation  of  (g^  =  0)  is  the  loss  evaluated  at  the  apex  of  the  Mth  (2nd) 

screen,  an  exact  value  of-201og|o(l/(M  =  2))  =  6.0206  dB  which  is  correctly  suggestive  of  the 

well  know  result  for  a  single  (M-1)  knife  edge  which  is  also  plotted  periodically  with  values 
shown  as  solid  black  circles.  The  exact  confluence  of  multi-screen  diffraction  theory  with  known 
results  for  a  single  knife  edge  inspires  the  inclination  to  attempt  a  general  diffraction  scheme 
involving  both  formalisms.  Direct  comparison  with  the  conventional  Walfisch-Ikegami  result  is 
shown  in  Figure  9.1  where  the  loss  at  900  MHz  is  plotted  as  a  function  of  scaled  base  elevation 
(normal)  deviation  from  building  alignment.  The  exact  result  using  the  equations  of  Boesma  are 
shown  again  as  solid  circles.  The  numerical  simulation  and  conventional  Walfisch-Ikegami  are 
observed  to  agree  well  with  the  exaet  result  for  positive  values  of  the  argument  but  the 
conventional  result  diverges  for  base  elevations  below  roof-top  alignment.  The  linear  behavior  of 
the  loss  below  roof  top  level  apparent  in  the  (COST  231)  Walfsich-Ikegami  model  (with  a  slope 
of  .8  dB/  meter)  is  characteristic  of  the  model’s  MSD  diffraction  loss  which  is  of  the  specific 
form. 


+  K  log  1 0{d)  +  log  1 0( V J  -  9  log  1  Q{b) 

where  for  base  station  elevations  above  roof  top  level,  k_a  is  given  by 

^,=54-.8(A/i,) 


(7A) 


(7B) 
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XS  LOSS  (DB)  vs  GC 


GC=DY/SQRT(wL  B) 

Figure  9 

MULTISCREEN  DIFFRACTlONiEXACT.CONVENTIONAL  AND  PROPOSED 


Figure  9.1 
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4.1. 1.2.  Extension  of  Urban  Vertical  Plane  Propagation  to  Non-Uniform  Building 
Heights  and  Inter-Building  Spacing 

The  proposed  model  also  considers  urban  terrain  profiles  more  complex  than  the  simple  two  ray 
paths  addressed  by  the  Walfisch-Ikegami  model.  Whereas  the  Walfisch-Ikegami  model  is 
restricted  to  urban  terrain  consisting  of  at  most  two  ray  paths  including  an  unobstructed  ray 
associated  with  multi  screen  diffraction  and  a  final  ray  associated  with  roof-to-street  diffraction, 
the  proposed  model  considers  more  general  urban  terrain  as  shown  in  Figure  10  which  may 
include  localized  segments  conforming  to  the  Walsich-Ikegami  assumptions,  but  may  also 
include  other  portions  which  are  not  aligned  with  a  segment  of  local  uniformity.  In  the  proposed 
model  the  partial  attenuations  attributable  to  such  portions  of  the  terrain  are  computed  according 
to  multiple  (optically  independent)  knife  edge  diffraction.  The  idealized  multi-screen  diffraction 
(MSD)  portions  are  indicated  in  Figure  10  with  white  arcs  spanning  the  angle  between  effective 
source  and  receiver  horizon.  The  independent  knife  edge  diffraction  centers  are  indicated  by 
black  arcs  centered  at  the  screen  which  subtend  the  angle  between  the  localized  source  and  its 
effective  receiver  horizon.  The  urban  terrain  of  Rosslyn,  Virginia  rarely  offers  profiles  in 
conformity  with  the  idealized  Walfishch-lkegami  assumptions  and  is  exemplified  by  the  terrain 
profile  shown  in  Figure  1 1. 


Figure  10 
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Figure  11 

4.1. 1.3.  Alternative  Formulation  of  Roof-to-Street  Diffraction  with  Recent 
Heuristic  Modifications  of  LITD 

A  final  perceived  deficiency  of  the  Walfisch-Ikegami  model  is  its  use  of  flawed  diffraction 
coefficients  from  a  dated  formulation  of  the  geometrical  theory  of  diffraction  which  have  since 
been  corrected.  In  its  modeling  of  roof-to  street  diffraction,  the  Walfisch-lkegami  model  uses  a 
diffraction  coefficient  for  perfectly  absorbing  wedges  and  screens  derived  by  Felsen  [2]  from 
geometrical  optics  and  therefore  incorrect  for  configurations  approaching  an  incidence  shadow 
boundary.  The  UTD  diffraction  coefficients  of  Kouyoumjian  and  Pathak  [4]  are  properly 
behaved  in  the  shadow  boundary  and  admit  to  an  exact  solution  for  screen  diffraction  (exterior 
wedge  angle  =  2n ).  In  addition  Holm  has  shown  that  his  heuristic  UTD  modeling  of  wedge 
diffraction  for  imperfectly  conducting  materials  reduces  to  a  single  partial  diffraction  coefficient 
in  the  perfectly  absorbing  limit  obtained  from  the  following  expression  for  =  0 , 

D  =  D,  +  +{R^R^  )Dj  +  R„D^  (8) 

Such  a  formulation  of  roof-to-street  diffraction,  generalized  to  include  reflection  from  the  walls 
of  adjacent  buildings,  is  employed  in  the  proposed  model. 
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4.1.2.  Urban  Terrain  Geometrical  and  Other  Considerations  for 
Implementation  of  Proposed  Methods 

4.I.2.I.  Extraction  of  Terrain  Profile 

Provided  with  traditional  data  files  for  urban  terrain,  faces  of  constant  elevation  are  first  detected 
and  interpreted  as  floor  levels  or  rooftops  of  buildings.  The  surface  dimensions  of  the  vertices 
and  their  constant  elevation  are  stored  and  a  building  identification  index  is  assigned  to  the 
coordinates  so  stored  for  immediate  retrieval.  The  information  collected  can  be  used  to  create  an 
urban  building  footprint  map  as  appears  in  Figure  12.  Given  a  point  of  initiation  on  the  urban 
terrain  for  the  profile  of  interest  and  either  a  receiver  location  or  intended  orientation  of  the 
vertical  plane  relative  to  the  origin,  the  intersection  of  the  vertical  plane  with  the  buildings  lying 
in  its  path  may  be  computed  to  create  an  array  of  distances  from  the  transmitter  in  the  one 
dimensional  range  coordinate  and  the  vertical  elevation  of  the  points  of  intersection  along  the 
profile  at  each  intersection.  The  intersections  in  the  footprint  map  appear  in  the  Figure  12  as 
solid  (black)  circles  and  the  corresponding  urban  profile  in  the  resulting  vertical  plane  is  shown 
in  Figure  13.  For  the  purposes  of  examining  an  urban  terrain  profile  exhibiting  a  more  general 
composition  including  localized  structures  consistent  with  multi-screen  diffraction  as  well  as 
optically  isolated  knife  edges,  the  urban  terrain  profile  shown  in  Figure  14  will  be  considered  in 
the  discussion  which  follows. 

APICAL  FACE  TANGENTIAL  PROJECTION 


Figure  12 
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VERTICAL  PLANE  URBAN  PRORLE 


Figure  13 

MODEL  URBAN  TERRAIN  PRORLE 


D(M) 

Figure  14 
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4.I.2.2.  Formalism  for  Geometrical  Deflnition  and  Detection  of  Buildings 

Since  as  suggested  by  Bertoni  [7]  [8],  buildings  are  to  be  approximated  by  single  diffracting 
screens,  it  becomes  necessary  to  establish  a  formalism  for  geometrically  defining  a  building. 
Aware  of  the  definition  of  a  diffraction  centers  as  any  local  convexity  introduced  by  physical 
structures  to  incident  ray  paths,  the  definition  of  a  building  adopted  for  this  study  is  any  region  of 
monotonic  local  curvature  in  which  the  inclination  of  ray  paths  is  reduced  by  each  successive 
component  vertex.  The  building  identification  sequence  resulting  form  this  definition  for  the 
urban  terrain  shown  in  the  Figure  14  as  red  numerals  corresponding  to  building  count 
assignments.  A  corresponding  sequence  of  index  arrays  are  created  for  identifying  the  terrain 
index  associated  with  a  particular  building  assignment  and  an  inverse  array  for  recovering  a 
possible  building  assignment  associated  with  a  given  terrain  index.  The  accumulation  of  the 
building  identification  sequence  permits  the  construction  of  individual  building  parameter  arrays 
which  contain  the  coordinates  (range  and  heights)  of  the  vertices  understood  to  comprise  a  given 
building.  This  information  is  vital  to  establishing  building  statistical  parameters  such  as  average 
width  and  periodic  spacing,  pending  identification  of  the  relevant  multi-screen  diffraction 
configurations.  Anticipating  the  recommended  constraint  that  each  building  provide  but  a  single 
absorbing  screen  in  the  effective  urban  diffraction  model,  the  urban  terrain  profile  is  next 
reduced  to  a  sequence  of  half  screens  of  negligible  thickness,  resulting  in  the  essential  urban 
profile  shown  in  Figure  15. 


GEODESIC  PATHS 
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Figure  15 
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4.1.2.3.  Extraction  of  Indexed  Geodesic  Path  to  all  Points  in  Terrain  Proflle 

The  generalized  Fermat  Principle  for  diffraction  states  that  a  diffracted  ray  path  will  proceed 
from  the  transmitter  to  the  receiver  along  the  shortest  possible  (geodesic)  path.  In  anticipation  of 
the  most  general  preference  to  evaluate  path  losses  at  an  arbitrary  intermediate  range  of  the 
profile  and  at  an  arbitrary  height,  and  in  view  of  Fermat’s  principle  for  diffraction,  it  becomes  of 
interest  to  systematize  a  deseription  of  the  points  encountered  in  a  geodesic  excursion  from  a 
fixed  transmitter  (understood  to  oeeupy  a  specified  elevation  at  the  beginning  of  the  profile)  to 
any  intermediate  point  in  the  profile.  An  intermediate  profile  locus  is  understood  to  oeeupy  a 
nominal  indexed  vertex  of  the  profile  or  is  defined  as  a  linear  interpolation  between  any  two  such 
vertices.  Appealing  to  the  definition  of  geodesic  paths  in  a  plane  as  defining  a  sequence  of  line 
segments  proceeding  through  a  sequence  of  monotonically  deereastng  slopes,  an  algorithm  was 
written  whieh  identifies  all  such  paths,  and  records  the  path  index  to  all  points  in  the  profile  and 
all  constituent  points  encountered  in  a  given  path.  A  schematic  presentation  of  the  various  paths 
to  all  points  in  the  reduced  profile  is  shown  in  Figure  15. 

4.1.2.4.  Accumulation  of  Geodesic  Paths  and  Horizon  Arrays  for  all  Points  in 
Path 

Provided  with  the  indexed  paths  and  path  contents  to  an  arbitrary  intermediate  terrain  location 
for  which  the  path  loss  is  desired,  it  becomes  necessary  to  examine  the  properties  of  the 
associated  path  and  to  identify  and  record  diffraction  centers  encountered  along  the  path  with 
proper  regard  to  elassification  as  either  a  component  of  a  multi  screen  diffraction  construct  or  as 
an  element  of  an  optically  independent  knife  edge  sequence.  Since  the  path  to  a  given  point  is 
geodesic  and  the  geodesic  paths  are  as  depicted  in  Figure  15,  the  intermediate  points  of  the  path 
to  a  specified  point  eaeh  define  a  center  of  diffraction  along  the  path  to  the  static  receiver.  As  an 
example  the  path  to  the  final  point  of  the  profile  is  shown  in  Figure  16  where  the  geodesic  path  is 
drawn  in  red  with  the  intermediate  diffraction  centers  indicated  by  partial  arcs  centered  at  each 
diffraction  center  and  spanning  the  angle  between  the  corresponding  effective  transmitter  and 
receiver. 
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LOCAUZED  MULTI-SCREEN  DIFFRACTION 


Figure  16 


4.I.2.5.  Identifications  of  MSD  Diffraction  Centers 

The  optical  environment  of  each  ray  of  the  geodesic  path  is  next  inspected  to  determine  whether 
any  intervening  (sub-geodesic)  structures  penetrate  the  first  Fresnel  zone  between  the  endpoints 
of  the  geodesic  ray  considered.  Such  sub-geodesic  intrusions  are  readily  apparent  in  Figure  1 6  as 
knife  edges  penetrating  the  relevant  first  Fresnel  zones  which  are  drawn  for  each  segment  of  the 
geodesic  path.  Thus  the  first  Fresnel  zone  between  points  (0)  and  (24)  in  the  figure  is  penetrated 
by  profile  points  (9)  and  (12).  The  sub-geodesic  points  of  intrusions  and  the  final  profile  point  of 
the  segment  are  interpreted  as  defining  a  multi-screen  diffraction  configuration  in  which  the 
corresponding  building  dimensions  are  used  to  define  the  average  spacing  and  building  width  as 
required  to  compute  the  relevant  MSD  contribution  to  the  excess  path  loss.  Designating  the  final 
point  of  the  path  as  a  center  of  rotation,  the  average  orientation  of  roof  tops  is  defined  as  the 
average  angle  from  this  appointed  origin  to  each  of  the  optically  relevant  buildings.  The  multi¬ 
screen  diffraction  loss  is  also  a  function  of  the  normal  coordinate  above  the  average  roof  top 
alignment,  the  number  of  knife  edges  found  to  penetrate  the  first  Fresnel  zone  and  their  average 
spacing.  The  MSD  constructs  are  identified  by  three  terrain  indices  which  are  the  effective 
transmitter  (first  point  of  a  geodesic  ray),  the  first  knife  edge  found  to  be  of  optical  relevance  and 
the  endpoint  of  the  geodesic  ray.  Three  such  MSD  diffraction  centers  were  identified  in  the 
example  urban  terrain  and  are  shown  in  Figure  16  as  red  arcs  centered  on  the  second  of  these 
identifying  terrain  features  and  spanning  the  angle  between  the  relevant  ray  path  endpoints. 
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4.I.2.6.  Consolidation  of  Previously  Independent  MSD  Constructs 

Whereas  the  methodology  described  above  activates  sub-geodesic  intrusions  between  principle 
vertices  of  a  given  geodesic  path,  an  auxiliary  procedure  involves  inspection  of  adjacent  geodesic 
vertices  for  sufficient  optical  alignment  to  justify  incorporation  of  adjacent  elements  into  a  single 
multi  screen  diffraction  center.  The  common  vertex  of  two  adjacent  but  independent  diffraction 
centers  is  then  redefined  as  a  diffracting  element  of  a  single  multi-screen  configuration 
composed  of  these  originally  independent  elements.  Thus  the  first  two  optical  diffraction 
constmcts  shown  in  the  top  panel  of  Figure  17,  initially  modeled  as  (1)  a  single  knife  edge 
involving  the  transmitter  and  first  two  dif&action  screens  followed  by  (2)  a  19  screen  multi¬ 
screen  diffraction  configuration,  are  consolidated  to  create  the  single  MSD  shown  in  the  second 
panel  of  Figure  17. 
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Figure  17 


4.1.3.  Loss  Evaluation  Process 

4.I.3.I.  Evaluation  of  Path  Loss  at  Arbitrary  Range,  and  HI  and  H2 

Whereas  the  earlier  characterization  of  the  geodesic  path  considered  the  terminal  range  of  the 
urban  terrain,  the  quantitative  inclusion  of  comparable  data  for  reconstructing  the  geodesic  path 
to  an  intermediate  range  receiver  makes  possible  the  evaluation  of  excess  path  loss  at  any  such 
point.  Since  the  principle  diffraction  centers  are  the  vertices  of  geodesic  paths  between  antennas, 
provision  of  the  geodesic  path  provides  geometrical  information  essential  to  the  calculation  of 
path  losses  involving  diffractive  mechanisms. 
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To  illustrate  the  mechanism,  an  elevated  receiver  at  an  intermediate  range,  as  labeled  in  Figure 
18,  is  analyzed  by  first  recalling  the  (receiver)  horizon  of  its  most  immediate  terrain  profile 
relation  in  the  direction  of  the  transmitter.  If  the  slope  from  the  most  immediate  horizon  exceeds 
that  of  the  geodesic  approach  to  the  horizon  from  the  transmitter,  the  geodesic  path  to  the 
horizon  is  examined  in  reverse  order  to  establish  the  horizon  of  the  specified  field  point.  Since 
the  geodesic  paths  to  all  points  in  the  profile  are  know,  identification  of  the  receiver  horizon 
among  the  profile  points  equivalently  identifies  the  geodesic  path  to  the  specified  (ground  level) 
field  point. 

In  the  absence  of  back-scattering  the  geodesic  path  to  ground  level  provides  a  fiill  description  of 
the  optical  environment  of  a  receiver  of  any  height  above  ground  level  and  in  general  is  required 
to  enforce  continuity  of  excess  loss  as  a  function  of  elevation  at  the  specified  range.  Provision  of 
the  geodesic  path  to  the  specified  receiver  permits  itemization  of  the  relevant  diffraction  centers 
which  are  examined  for  prior  classification  as  MSD  or  independent  knife  edge  diffraction 
centers.  If  an  intermediate  diffraction  center  is  identified  as  the  central  vertex  of  a  multi-screen 
diffraction  construct,  the  (multiple)  knife  edge  count  must  be  altered  appropriately  from  the 
initial  values  established  on  the  occasion  of  its  initial  detection.  Thus  the  knife  edge  multiplicity 
count  of  the  MSD  considered  MSD  (2,3,4)  in  Figure  18  must  be  altered  from  its  original  value  of 
four  to  the  two  intervening  knife  edges  between  the  central  vertex  at  650  meters  and  the  specified 
receiver.  Alteration  of  the  geodesic  vertex  assignments  for  the  intermediate  receiver  may  also 
alter  the  third  vertex  assignment  of  pre-existing  single  knife  edges  as  required  to  compute  altered 
maximum  attenuation  bounds  appropriate  for  intermediate  shadows.  The  final  itemization  of 
modified  diffraction  centers  appropriate  for  calculation  of  the  excess  loss  at  the  (elevated) 
intermediate  field  point  shown  in  Figure  18  includes  the  first  MSD  diffraction  center 
MSD(0,1,2),  a  single  knife  edge  centered  at  point  (2)  and  redirected  to  the  intermediate  receiver 
horizon  and  a  second  multi-screen  diffraction  center  MSD(2,3,RXH)  reconfigured  for 
termination  at  the  modified  receiver  horizon  with  a  reduced  constituent  knife  edge  count  of  two. 

Before  detailing  the  numerical  accumulation  of  the  various  partial  excess  losses,  it  is  appropriate 
to  recall  the  approximations  of  multiple  knife  edge  diffraction  originally  due  to  Epstein  and 
Peterson  [9]  which  approximates  the  total  multi-screen  diffraction  loss  as  a  sum  of  single  knife 
edge  attenuation  subject  to  an  upper  bound  imposed  by  any  subsequent  diffraction  screens.  A 
typical  urban  terrain  profile,  after  reduction  to  its  component  screens  and  consisting  only  of 
optically  independent  knife  edges,  is  shown  in  Figure  19  with  schematic  arcs  indicating  the 
angular  extent  of  each  knife  edge.  At  relatively  remote  elevations,  the  second  knife  edge 
introduces  a  negligible  attenuation  as  suggested  by  the  black  vertical  line  segment  in  the  figure. 
The  cascaded  implementation  of  successive  knife  edges  requires  that  each  component  be 
deactivated  at  a  specific  maximum  attenuation  as  the  receiver  approaches  an  elevation  where  the 
succeeding  knife  edge  begins  to  exert  an  attenuation.  This  critical  configuration  corresponds  to 
the  elevation  at  which  the  succeeding  knife  edge  first  penetrates  the  first  Fresnel  zone  between 
its  effective  source  and  the  variable  receiver  height.  Thus  the  attenuation  of  a  given  knife  edge  is 
set  equal  to  a  constant  at  the  critical  elevation  corresponding  to  the  onset  of  attenuation  of  its 
geodesic  successor.  Considering  the  succession  between  first  and  second  knife  edges  in  the 
Figure  19,  the  locus  of  points  beyond  the  second  knife  edge  associated  with  critical  penetration 
of  the  relevant  Fresnel  zone  (between  the  first  knife  edge  and  variable  receiver)  is  shown  as  a 
blue  dashed  curve  extending  from  the  second  knife  edge. 
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Figure  18 
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4.I.3.2.  Final  Accumulation  of  Excess  Loss 

Returning  to  the  final  accumulation  of  excess  loss  appropriate  for  the  receiver  located  at  the 
intermediate  range  of  Figure  19,  the  succession  of  independent  excess  losses  due  to  each 
diffraction  center  is  made  apparent  by  considering  an  excursion  in  elevation  initiated  well  above 
the  line  of  sight  beyond  any  diffractive  influences.  As  the  receiver  elevation  approaches  the 
greater  threshold  of  the  first  multi-screen  diffraction  center,  a  component  of  loss  due  to  the  MSD 
is  introduced  which  increases  to  the  maximum  value  obtained  at  the  roof  top  level  of  the  final 
knife  edge  of  the  MSD  construct.  The  excess  loss  continues  to  increase  as  the  second  diffraction 
center,  SKE  (0,2,RXH)  introduces  a  diffractive  loss  between  its  initiation  threshold  and  that  of 
the  final  knife  edge  centered  at  the  receiver  horizon.  In  the  final  segment  of  the  excursion  and 
descent  to  ground  level,  the  diffractive  influence  of  the  final  diffraction  center  SKE  (2,  RXH, 
RX)  is  added  to  the  excess  loss  resulting  in  a  total  excess  loss  at  the  receiver  of  52.3  dB.  The 
various  components  of  the  total  excess  loss  and  their  sum  have  been  plotted  as  a  function  of  the 
structural  receiver  antenna  height  in  Figure  20. 
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EXCESS  LOSS  AT  INTERMEDIATE  RANGE 


Figure  20 


4.2.  Fast  Urban  Predictions 

4.2.1.  Initial  Effort 

For  the  initial  effort,  it  was  deemed  prudent  to  follow  the  methods  used  by  MITRE  in  [10].  This 
effort  proved  to  be  a  starting  point  for  a  fast  urban  path  loss  estimator.  The  optimized  ray  engine 
(ORE)  was  used  in  conjunction  with  a  building  coimting  algorithm  and  a  path  loss  attenuation 
formula.  This  algorithm  computes  path  loss  based  on  the  amount  of  free  space  path  loss 
combined  with  path  loss  due  to  urban  settings.  Computing  reflections,  transmissions  and 
diffractions  aren't  necessary  because  the  building  depth  accounts  for  their  contribution.  The  first 
implementation  step  used  ORE  to  count  the  number  of  facets  between  the  transmitting  and 
receiving  antennas. 
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Figure  21:  Colors  showing  number  of  walls  passed  through 

From  there,  a  simple  algorithm  calculated  the  number  of  buildings  between  the  transmitting  and 
receiving  antenna. 
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Figure  22:  Colors  showing  number  of  buildings  passed  through 

Now,  knowing  the  building  depth  and  the  free  space  depth  formulas  from  the  [10]  could  be  used 
to  calculate  path  loss  attenuation.  Further  work  included  automatic  tuning  of  the  model  based  on 
building  material  parameters  and  frequency  as  well  as  adding  additional  propagation 
mechanisms. 
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Figure  23:  Left  Wireless  InSite  Urban  Canyon  2D,  Right  new  model 


Figure  24:  Left  Urban  Canyon  FDTD,  Right  new  model 

Further  Remcom  compared  the  results  of  this  model  with  both  measured  data  and  Wireless  InSite 
predictions.  These  are  shown  in  Figure  25,  Figure  26,  Figure  27  and  Figure  28. 
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Figure  25:  Comparisons  of  Path  Gain  on  Albert  Street,  Ottawa  Canada 
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Figure  26:  Comparisons  of  Path  Gain  on  Bank  Street,  Ottawa  Canada 


Contract  Number  W3 1P4Q-04-C-R245 


35 


14  August  2006 


Figure  27:  Comparisons  of  Path  Gain  on  Laurier  Street,  Ottawa  Canada 


Figure  28:  Comparisons  of  Path  Gain  on  Queen  Street,  Ottawa  Canada 
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Due  to  the  semi-empirical  nature  of  the  fast  prediction  methods  validation  and  tuning  of  the 
models  are  of  great  importance.  Initial  testing  demonstrated  in  the  above  figures  revealed  the 
need  to  incorporate  higher  order  propagation  phenomena  such  as  diffi'action  and  reflection  to 
account  for  discrepancies  in  long  range  paths  that  are  dominated  by  such  phenomena.  How  these 
effects  and  other  geometric  properties  of  the  urban  environment  such  as  material  definition, 
frequency,  transmitter  and  receiver  height,  building  geometry,  and  antenna  pattern  must  all  be 
factored  into  the  model  as  appropriate  weighting  terms  on  loss  functions  has  been  the  focus  of 
current  research.  A  series  of  tests  was  designed  using  REMCOM’s  Wireless  InSite  tool  to  study 
the  effects  of  these  properties  using  a  GTD/UTD  based  urban  canyon  model  and  an  Urban 
Canyon  FDTD  model.  The  results  of  these  tests  were  then  used  as  a  baseline  along  with 
available  measured  results  to  tune  the  weighting  parameters  of  the  fast  prediction  model.  Figure 
29  and  Figure  30  some  of  the  higher  order  effects  that  can  be  determined  by  our  fast  model. 


Figure  29  Multiple  over  the  roof  top  diffraction  found  by  fast  model 
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Figure  30  Multiple  effects  (horizontal  plane  diffraction,  vertical  over  the  roof  top  diffraction,  and 

single  reflections) 


4.2.2.  Technical  Details 

The  ORE  engine  makes  use  of  spatial  subdivisioning  techniques  and  parametric  traversal  to 
compute  a  path  through  the  urban  environment.  This  path  is  composed  of  segments  that  are 
characterized  by  the  media  through  which  they  are  traversing.  Each  segment  that  occurs  within 
an  obstruction  is  attenuated  using  a  statistical  method  based  on  its  geometric  characteristics. 
These  segments  are  combined  to  provide  a  final  path  loss  value  on  the  order  of  1000/second.  If 
more  time  is  available  the  engine  will  be  capable  of  finding  critical  first  order  interactions  and 
including  these  into  the  approximation  for  increased  accuracy. 

4.2.3.  MITRE  Model 

The  MITRE  model,  [10],  was  selected  to  be  the  first  model  to  test  to  meet  the  run  time  criteria. 
There  are  two  tuning  parameters  available  in  this  model.  The  attenuation  due  to  buildings  is 
described  by  the  formula  below. 


Lb  =  K,  *(0.6/+K2*Db 

Db  is  the  amount  of  building  depth  in  meters;  f  is  Frequency  in  GHz  and  the  two  K  factors  are 
used  to  tune  the  building  attenuation. 
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Given  measured  or  calculated  results,  K  factors  can  be  tuned  using  a  best  fit  line  solution.  Two 
approaches  were  available  to  alter  these  factors  in  each  given  scenario.  If  the  factors  were 
consistent  throughout  the  given  urban  geometry,  a  high  fidelity  model  could  run  as  a 
preprocessing  phase.  Measured  results  in  the  cities  of  Ottawa,  Bern,  and  Rosslyn  were  used  in 
conjunction  with  Wireless  InSite  UTD  calculated  results  to  determine  the  consistency  of  the 
behavior  in  a  given  urban  scenario.  Below  are  graphs  showing  K1  and  K2  vs  Transmitter  height 
and  frequency  in  order  to  determine  if  they  are  constant  on  a  per  city  basis. 


K1  vs  Tx  Height 
In  Rosslyn,  VA 
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K1  vs  Frequency 
in  Rosslyn,  VA 


K2  vs  Frequency 


Tx1  K2 
Tx2  K2  ^ 
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K1  vs  Receiver  Height 
in  Ottawa,  CA 
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K2  vs  Receiver  Height 
in  Ottawa,  CA 
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K1  vs  Frequency 
in  Ottawa,  CA 


Once  the  analysis  was  completed  it  was  clear  that  the  factors  were  not  consistent  throughout  the 
city  and  no  correlation  could  be  determined  based  on  relative  transmitter  or  receiver  location. 

Since  the  high  fidelity  tuning  failed,  our  second  approach  was  to  tune  the  MITRE  model  using 
the  empirical  COST-HATA  model.  The  theory  being  that  COST-HATA  was  relatively  good  for 
certain  scenarios,  and  could  be  more  deterministic  knowing  building  depth.  Below  are  results  for 
three  eities  using  the  MITRE  model  tuned  with  the  empirical  model. 
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Ottawa 

K1  =37.9877 
K2  =  .0393 


RMS  Comparisons  to  WI 
LIPPS;6.6678 
CostHata:9.0588 


RMS  Comparisons  to  WI 
UPPS:  10.78 
CostHata:!  1.566 


RMS  Comparisons  to  WI 
UPPS:  10.7795 
CostHata:  10.2 139 
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RMS  Comparisons  to  WI 
UPPS:6.4895 
CostHata:7.435 


Rosslyn 

K1  =  39.66 
K2  =  .0603 
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Bern 

K1  =  80.79 
K2  =  .0925 
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As  seen  from  the  graphs,  the  results  are  better  than  COST-HATA,  but  mostly  only  in  line-of- 
sight  regions.  The  conelusion  from  these  tests  is  that  tuning  the  MITRE  model  on-the-fly  does 
not  have  the  accuracy  desired.  Another  drawback  to  using  an  empirical  model  to  tune  the 
MITRE  model  is  that  empirical  models  are  only  valid  in  limited  environments. 

From  these  two  studies  it  became  apparent  that  the  MITRE  model  was  not  a  viable  solution  for  a 
fast  accurate  deterministic  model.  Efforts  shifted  towards  researching  other  possible  models. 
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4.2.4.  Other  models 

The  C-style  API  was  used  to  integrate  the  newly  developed  fast  urban  models  into  network 
modeling  tools  such  as  OPNET  and  QualNet.  Fiuther  an  API  manual  was  completed.  In  order  to 
facilitate  using  the  API,  a  sample  program  was  written  as  an  example  for  first  time  users.  The 
API  includes  the  following  empirical  propagation  models: 

HATA 

Cost-HATA 

Freespace 

OPAR 

Extended  Walfisch-Ikegami  (WI)  model. 

Further  to  make  it  easier  to  use  some  of  the  models  algorithms  have  been  implemented  to 
determine  city  type  for  HATA,  Cost-HATA  and  WI  depending  on  building  heights  and  building 
densities. 

As  part  of  an  effort  to  implement  the  vertical  pane  urban  path  method  model,  the  new  ray  tracer 
was  enhanced  to  quickly  query  the  status  of  a  point,  whether  it  is  outside,  inside  a  building  or 
inside  the  terrain.  This  method  uses  geometry  pre-processing,  edge  caching,  a  method  for  fast 
vertical  and  horizontal  convex  hull  determination.  As  the  dominant  path  method  is  refined  and 
the  algorithm  understood  better,  this  method  shall  be  included  in  the  API. 

4.3.  Building  Simplification 

4.3.1.  Two  Dimensional 

4.3.1. 1.  Simplification  Technical  Details 

A  common  problem  for  urban  propagation  is  generating  appropriate  representation  of  the 
buildings.  Too  many  faces  on  the  buildings  cause  anomalous  propagation  results  as  well  as 
extensive  computation  time.  This  will  be  especially  troublesome  when  urban  models  are 
hybridized.  Thus,  an  easy  to  use,  transparent  to  the  user,  and  computationally  efficient  method  is 
required  to  simplify  buildings.  Remcom  has  investigated  many  methods  fi'om  computer  graphics 
to  solve  this  problem.  Recently  Remcom  tested  an  improved  method  for  urban  canyon 
propagation.  This  method  successfully  reduced  the  city  shown  in  Figure  31  to  one  less  than  half 
the  original  size  (Figure  32)  and  finally  to  a  final  rendition  (Figure  33)  that  was  10%  of  the 
original  size. 
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Figure  31:  DuPont  Circle:  Original  (30,000+  vertices) 
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Figure  32:  DuPont  Circle:  first  level  of  simplification  (14,430  vertices) 
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Figure  33:  DuPont  Circle:  flnal  simpliflcation  (3,135  vertices) 


The  simplified  city  models  were  then  used  to  compute  propagation  loss  (Figure  35)  and  the 
results  were  compared  to  the  results  found  for  a  full-resolution  city  (Figure  34). 


Figure  34:  DuPont  Circle  propagation  coverage  example  with  full  building  resolution 

In  many  areas  the  propagation  loss  is  comparable.  In  fact  the  results  seem  to  differ  most  in  the 
inner  courtyards.  This  is  largely  caused  by  the  unioning  groups  of  buildings  into  one. 
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Figure  35:  Example  of  the  above  test  with  the  simplified  city. 

To  better  quantify  the  differenees  the  loss  was  computed  along  the  two  routes  shown  in  Figure 
36,  Figure  37,  and  Figure  38. 


Figure  36:  Two  example  routes  used  to  compare  the  results 
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Figure  37:  The  loss  down  the  diagonal  street  for  the  original  city  (blue)  and  the  simplified  city  (red) 


Figure  38:  The  loss  along  the  cross  street  for  the  original  city  (blue)  and  the  simplified  city  (red) 

These  methods  to  simply  building  footprints  in  support  of  urban  canyon  modeling  were 
integrated  into  Wireless  InSite®.  The  previous  method  applied  sophisticated  methods  but  was 
best  applied  to  buildings  of  similar  complexity.  That  is,  often  the  tolerances  for  simplification  for 
a  complex  building  were  too  stringent  for  a  simple  building.  Thus,  global  metric  and 
corresponding  algorithms  are  being  developed  to  allow  the  user  to  use  a  single  metric  to  achieve 
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varying  degrees  of  building  simplification  over  a  large  portion  of  a  city.  This  involves  balancing 
the  reduction  of  simpler  shaped  buildings  with  those  that  are  more  complex.  The  Douglass- 
Pueker  polyline  simplification  algorithm  was  modified  to  be  performed  on  polygons.  Utilization 
of  this  algorithm  provides  a  global  metric  that  is  based  on  local  variants.  In  our  case  a  tolerance 
was  constructed  for  simplification  of  a  vertex  within  the  building  footprint  using  a  percentage  of 
the  distance  from  the  vertex  to  the  footprint  centroid.  This  tolerance  acts  as  the  tolerance  found 
in  Douglass  and  Pucker’s  initial  work.  A  more  in  depth  sensitivity  study  shall  be  constructed  in 
order  to  validate  that  the  building  simplification  does  not  unduly  perturb  the  propagation 
predictions. 

4.3. 1.2.  Integration  in  Wireless  InSite 

The  simplification  library  was  utilized  within  Wireless  InSite  to  perform  building  simplification. 
The  current  simplification  process  opens  a  window  as  shown  below. 
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The  window  has  two  sliders;  one  for  controlling  the  amount  of  simplification  to  perform,  and 
another  that  controls  the  amount  of  buildings  to  remove  based  on  established  criteria.  The 
Complexity  slider  controls  the  tolerance  parameter  for  the  Douglass-Puecker  simplifier  that  is 
used.  Lowering  the  complexity  raises  the  amount  of  tolerance  that  is  allowed  and  therefore 
removes  more  points  from  the  final  product. 

The  Building  Retention  slider  will  control  the  minimum  size  that  a  building  must  have  in  order  to 
be  retained.  At  present,  retention  can  be  determined  by  height,  area  or  volume.  It  was  decided  to 
keep  all  of  these  options  as  they  make  sense  for  different  models.  If  you  were  to  remove 
buildings  based  on  area,  then  narrow  towers  would  be  removed  that  would  be  important. 
However,  if  you  removed  buildings  based  on  volume,  then  the  tower  would  be  more  likely  to  be 
retained  over  lower,  broader  buildings  that  might  not  be  as  relevant  to  the  calculation.  Buildings 


54 


Prepared  by  Remcom,  Inc. 


Final  Report 


with  differing  materials  can  be  combined,  but  is  provided  as  a  deactivated  option  as  materials 
play  an  important  part  in  the  final  calculation. 

Information  is  provided  about  the  simplification  process  that  shows  how  much  of  an  impact  the 
current  settings  have  had.  The  user  can  also  decide  to  show  the  transmitter  and  receiver  sets  of 
the  project  in  the  simplification  window.  This  is  useful  for  determining  if  the  simplification 
process  caused  a  building  to  envelope  any  or  all  points  of  a  set,  in  which  case  the  client  may 
decide  to  increase  the  complexity  in  order  to  prevent  this  occurrence. 

Here  is  an  example  of  the  simplification  performed  on  the  above  geometry: 


4.3.2.  Three  Dimensional  Simplification 

Our  efforts  on  building  simplification  focused  on  building  foot  prints  in  a  2  and  2.5D  scenario. 
Although  this  method  alone  has  shown  promising  results  with  urban  canon  FDTD  testing  it  was 
sufficient  to  attain  acceptable  results  and  significantly  reduce  run  times  when  using  full  3D  ray 
tracing.  The  problem  of  full  3D  simplification  is  very  difficult.  Many  areas  of  geometric 
modeling,  computer  graphics  and  finite  element  analysis  have  developed  innovative  algorithms 
for  3D  triangular  manifold  simplification.  These  algorithms  however  are  highly  constrained  on 
representations  of  geometry  to  fully  connected  well  formed  triangular  meshes  and  often  primary 
concern  is  with  visual  appearance  of  the  model  instead  of  the  electromagnetic  complexity. 
Couple  this  with  the  fact  that  urban  models  very  rarely  are  represented  within  these  constraints, 
and  to  modify  the  models  to  conform  is  an  equally  challenging  problem,  and  many  known 
methods  become  unusable.  Therefore,  research  in  this  area  was  focused  on  more  novel 
approaches  that  allow  flexibility  in  model  representation  while  reducing  runtime  requirements 
and  maintaining  accuracy  in  EM  simulation. 
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One  such  technique  investigated  was  introduced  by  Josilin,  this  technique  attempts  to  retrieve  the 
significant  facet  from  a  group  of  neighboring  facets  in  the  model  geometry.  This  significant 
facet  is  then  used  to  represent  the  entire  group.  One  of  the  most  attractive  aspects  of  the  method 
introduced  is  that  it  relies  on  spacing  parameters  that  can  be  modified  to  have  frequency 
dependence.  This  would  allow  custom  model  simplification  based  on  frequency  to  enhance  the 
accuracy  of  results.  Other  methods  investigated  include  3D  feature  extraction  and  removal. 

The  extraction  of  features  was  done  using  a  pseudo  Constructive  Solid  Geometry  (CSG) 
approach.  Once  features  are  extracted  they  can  be  stored  into  a  database  according  to 
electromagnetic  size.  In  this  scenario  a  Level  of  Detail  (LOD)  model  could  be  constructed  to 
only  include  the  significant  features  at  a  given  frequency.  The  difficulty  with  this  method  is  the 
extraction  of  features  which  is  very  computationally  intensive  and  not  well  defined  in  the 
literature. 


4.4.  Extended  Semi-Empiricism  for  (Complex)  Soil  Permittivity 

4.4.1.  Preface 

The  initial  approach  to  computing  the  soil  dielectrics  used  a  bi-cubic  spline  method.  However,  it 
was  recognized  that  the  bi-cubic  spline  methodology  developed  during  the  initial  phases  of  this 
effort,  while  invoking  the  advantages  of  established  two-dimensional  interpolation  techniques, 
sacrificed  the  full  generality  of  a  working  algorithm  with  explicit  dependence  upon  soil  solids 
gravimetric  fractions.  The  current  algorithm  features  a  refinement  of  an  earlier  intermediate 
process  which  retains  the  explicit  dependence  upon  soil  solids  composition  and  uniquely  permits 
a  recent  modification  to  include  a  composite  of  two  established  semi-empirical  algorithms 
intended  for  the  frequency  ranges  300-1300  MHz  and  1400  MHz  -18GHz.  The  current  report 
summarizes  the  essentials  of  the  final  algorithm  with  a  particular  emphasis  upon  the  adaptation 
of  semi-empirical  methods. 

4.4.2.  Introduction 

A  comprehensive  numerical  algorithm  capable  of  evaluating  complex  permittivity  of  soil  for 
arbitrary  frequency,  volumetric  water  content,  and  soil  solids  composition  is  constructed.  The 
semi-empirical  methods  of  Hallikainen  and  Dobson  [13]-[14](1400  MHz-18GHz)  and  of 
Peplinksi  [15]-[16]  (300  MHz- 1300  MHz)  are  extended  to  lower  frequencies  (.01  MHz  -50 
MHz)  by  numerical  modeling  of  measured  data  and  supplemented  by  other  interpolative 
algorithms  to  uniformly  span  the  intended  frequency  range  of  .01  MHz  to  18  GHz.  The  log  of  the 
permittivity  as  a  function  of  frequency  is  displayed  in  Figure  39  which  is  labeled  along  the 
frequency  axis  to  indicate  the  principal  subdivisions  within  which  the  dominant  data  sources  and 
methodologies  determine  the  corresponding  computational  logic.  As  indicated  in  the  figure,  the 
semi-empirical  methods  of  Hallikainen,  Dobson  and  Peplinksi  [13],  [14],  [15]  and  [16]  are 
utilized  in  the  frequency  range  v  >  300  MHz.  In  an  intermediate  frequency  range  1-50  MHz,  the 
measurements  of  Campbell  [1 1]-[12]  collected  at  (1,2,5,  10,20  and  50)MHz  are  used  to  generate 
cubic  polynomials  in  volumetric  moisture  at  each  frequency.  The  use  of  the  polynomials  to 
calculate  complex  permittivity  in  the  range  1-50  MHz  is  then  extended  to  the  semi-empirical 
regime  with  an  interpolative  transition  between  50  and  300  MHz.  Finally  the  algorithm  is 
extended  to  frequencies  as  low  as  .01  MHz  with  the  measurements  of  Sternberg  [17]  for  sand 
and  Montmorollinite  clay  which  are  each  reported  at  single  values  of  moisture  from  .001  MHz  to 
the  lower  limit  of  the  second  frequency  range  reported  here  as  1 .0  MHz. 
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4.4.3.  Empirical  Algorithm  for  Frequency  Range  1400-18000  MHz 

The  original  work  of  Hallikainen  et  al  [13]  and  of  Dobson  [14]  provided  empirical  methods 
consisting  of  polynomials  quadratic  is  volumetric  moisture  at  selected  frequencies  between  1400 
MHz  and  18  GHz.  The  coefficients  of  the  quadratic  polynomials  were  modeled  as  linear 
combinations  of  gravimetric  fractions  of  the  primary  soil  components  silt,  sand,  and  clay.  Also 
provided  was  a  semi-empirical  algorithm  which  retained  the  desired  dependence  upon  accessible 
physical  soils  properties  such  as  volumetric  moisture,  bulk  and  specific  densities,  and  (weight 
fractions)  of  the  primary  soil  components  of  silt,  sand  and  clay.  The  explicit  form  of  the 
polynomials  provided  by  Hallikainen  and  Dobson  [13]-[14]  (parametrically  dependent  upon 

frequency  v  and  soil  solids  gravimetric  fractions  and  ,  and  explicitly  dependent  on 

volumetric  moisture  =  -^ )  is  given  by 

_  2  _ 

£  {{v,w^i,w^,,wJ,vJ  =  Y,c„vJ'  ( 1 ) 

n=0 

This  expression  for  e’  is  accompanied  by  an  identical  form  for  the  imaginary  part  of  the  complex 
permittivity  given  by 

_  2  _ 

e\{v,  w^.  cl  "  (2) 

/i=0 
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Each  coefficient  of  a  specified  power  n,  is  given  by  a  linear  combination  of  the  weight 
fractions  of  the  primary  soil  solids  sand  and  clay  and  of  the  explicit  form 


(3) 


where  the  dimensionless  scaled  masses  are  given  by  the  ratio  of  the  respective  component  mass 
to  that  of  (dry)  soil  according  to  the  expressions 


W,: 
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The  coefficients  of  like  power  from  each  of  the  varieties  of  soil  may  be  understood  to  occupy  the 
elements  of  a  column  vector  which  may  be  expressed  as  the  product  of  a  weight  fraction  matrix 
W  with  rows  corresponding  to  soil  varieties  and  (three)  columns  corresponding  to  each  of  the 
three  primary  soil  components  which  are  silt,  sand  and  clay. 

|c>„=W|a>„  (4) 


Although  theoretically  more  elegant  semi-empirical  formulations  have  been  presented,  the 
original  empirical  polynomials  ,  which  expressed  real  and  imaginary  parts  of  the  complex 
permittivity  as  quadratics  in  volumetric  moisture  (with  coefficients  linearly  dependent  upon  the 
soil  solids  gravimetric  fractions)  may  retain  their  appeal  for  simplicity,  replication  of 
experimental  data  and  dependence  upon  readily  available  soil  physical  quantities  and  because  of 
the  absence  of  a  uniform  semi-  empirical  algorithm  for  the  entire  frequency  range  of  300MHz- 
18GHz  This  report  will  therefore  consider  the  modeling  provided  by  both  the  (empirical) 
polynomials  and  of  the  corresponding  semi-  empirical  equations.  As  reported  by  Peplinski  [15]- 
[16]  the  semi-empirical  method  may  underestimate  the  real  permittivity  at  1400  MHz;  in 
addition  it  has  been  observ'cd  locally  that  the  equation  for  the  imaginary  part  of  the  permittivity 
at  1400  MHz  may  provide  non-physical  (negative)  results  for  soils  with  a  dominant  sand 
composition  and  at  large  values  of  volumetric  moisture.  Consequently  the  empirical  model 
presented  here  is  a  modification  of  the  polynomials  reported  by  Hallikainen  and  Dobson  in 
which  the  coefficient  of  the  sand  weight  fraction  in  the  expression  for  the  quadratic  coefficient  of 
the  imaginary  component  of  sand  at  1400  MHz  (c,=-.313  in  Table  II  of  reference  [13])has  been 

divided  by  four.  The  plots  shown  in  Figure  40  and  Figure  41  display  the  behavior  of  a  modified 
set  of  polynomials  for  the  five  varieties  of  the  soil  examined  by  Hallikainen  and  Dobson.  [13]- 
[14].  Also  shown  as  representative  of  measured  data  (as  open  circles)  are  results  obtained  from 
the  original  polynomials  of  Hallikainen  and  Dobson.  The  agreement  is  seen  to  be  identical  with 
the  exception  of  the  results  shown  for  1400  MHz  for  which  the  reported  coefficients  were 
modified  as  explained  above.  As  it  is  the  objective  of  this  work  to  provide  methods  for 
estimating  complex  permittivity  for  soils  of  arbitrary  solids  texture  composition,  including  those 
outside  the  loamy  compositions  considered  by  Hallikainen  and  Dobson,  Figure  41  displays  the 
behavior  of  the  complex  permittivity  modeled  by  the  empirical  polynomials  for  the  three  primary 
soils  represented  by  silt,  sand  and  clay  and  selected  to  be  of  a  composition  suggested  by  the  of 
work  Campbell  [1 1]-[12].  An  alternative  perspective  of  the  empirical  modeling  in  the  frequency 
range  1400MHz  -  18GHz  is  shown  in  Figure  42  and  Figure  43  which  display  (for  the  five  soils  of 
the  Hallikainen-Dobson  study  and  for  primary  soils  silt,  sand  and  clay)  both  components  of  the 
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complex  permittivity  as  a  function  of  the  logarithm  of  frequency  at  selected  values  of  volumetric 
moisture.  The  original  polynomials  of  Dobson  and  Hallikainen  are  used  with  the  exception  of 
1400  MHz  where  the  sand  dependent  coefficients  of  the  (quadratic)  coefficient,  c,  was  divided 

by  a  factor  of  four.  The  imaginary  part  of  the  permittivity,  shown  in  red,  is  observed  to  exhibit 
the  characteristic  minimum  above  1400  MHz  as  associated  with  the  lower  frequency  effects  of 
the  free  water  relaxation  near  20  GHz.  Below  1400  MHz,  the  curves  are  extended  by 
interpolation  to  the  highest  frequency  studied  by  Campbell  [11]-[12]  at  50  MHz.  (As  will  be 
discussed  more  thoroughly  in  later  sections,  the  measured  data  of  Sternberg  for  frequencies 
below  1  MHz  are  shown  in  the  upper  right  and  lower  left  panels  of  Figure  43,  corresponding 
respectively  to  sand  and  clay.  Both  sets  of  data  are  observed  to  fall  on  the  moisture  curve  which 
was  selected  for  the  corresponding  component  considered:  .095  for  sand  and  .706  for  clay.)  In 
the  following  section,  the  results  of  an  alternative  semi-empirical  treatment  are  presented  for  the 
same  selection  of  soil  texture  and  physical  parameters  and  frequency  range  as  considered  in 
Figure  40,  Figure  41,  Figure  42  and  Figure  43. 
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Figure  40 
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Figure  42 
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Figure  43 


4.4.4.  Independent  Semi-Empirical  Algorithms  for  Frequency  Ranges  300- 
1300  MHz  and  1400-18000  MHz 


The  initial  work  of  Hallikainen  and  Dobson  [13]-[14]  presented  a  set  of  parameters  for  a  semi- 
empirical  formalism  based  on  a  modified  version  of  the  Cole-Cole  procedure  for  computing 
permittivity  for  aqueous  suspensions  and  moist  inorganic  solids  intended  for  the  frequency  range 
1400  MHz  -18GHz.  Peplinksi  [15]-[16]  later  adapted  the  semi-empirical  method  for  the 
frequency  range  300-1300  MHz.  Figure  44  is  shown  a  representative  plot,  chosen  as  Field  II 
from  the  work  of  Hallkainen  and  Dobson,  of  the  entire  semi-empirical  frequency  range  of  300 
MHz  to  18  GHz  in  which  each  of  these  formalisms  is  implemented  in  its  respective  assigned 
range.  A  distinctive  discontinuity  is  apparent  in  the  intermediate  frequency  regime  between  1300 
and  1400  MHz.  Peplinski’s  modification  of  the  semi-empirieal  methodology  involved  a  linear 
scaling  and  shift  of  the  results  obtained  using  the  earlier  model  and  provided  an  alternate 
(optimized)  form  for  the  effective  conductivity. 
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Figure  44 


4.4.5.  (Unified)  Semi-Empirical  Algorithm  for  Frequency  Range  300- 
14000  MHz 

The  discontinuity  between  1300  and  1400  can  be  diminished  by  introducing  Peplinksi’s  modified 
forms  over  a  broader  range  (in  the  approach  to  the  transition  frequency  from  frequencies 
exceeding  1300  MHz)  and  with  switching  functions  to  introduce  the  changes  gradually  through  a 
transitional  continuum.  (Whereas  the  original  semi-empirical  algorithm  was  intended  for 
frequencies  exceeding  1400  MHz,  Peplinksi  extended  the  methodology  with  an  alternate  set  of 
parametric  equations  for  frequencies  in  the  range  300- 130a  MHz  .  In  the  present  work,  a 
composite  of  the  methods  has  been  used  for  the  full  frequency  range  of  300  MHz-18GHz).  The 
results  for  such  a  hybrid  semi-empirical  approach  appear  in  Figure  45  where  the  discontinuity  is 
essentially  undetectable  and  the  curves  smoothly  approach  the  two  bounding  limiting  forms 
appropriate  for  frequencies  above  or  below  1300  MHz.  Figure  46  presents  a  larger  frequency 
range  perspective  of  these  same  results  showing  the  interface  between  the  semi-empirical 
methods  above  300  MHz  and  the  empirical  polynomial  modeling  below  50  MHz.  Shown  in 
Figure  47  and  Figure  48  are  the  results  for  permittivity  as  a  function  of  volumetric  moisture 
using  the  semi-empirical  methodology  for  exactly  the  same  frequency  ranges  and  physical 
conditions  as  considered  in  the  plots  of  the  empirical  results  shown  in  Figure  40  and  Figure  41. 
Whereas  the  agreement  is  less  striking  than  that  obtained  with  empirical  polynomials,  the  semi- 
empirical  method  retains  the  appeal  of  its  more  theoretical  motivation  and  structure  and  its 
adaptability  to  similar  methods  at  lower  frequencies.  In  addition  the  data  points  shown  as  open 
circles  do  not  correspond  to  exact  measurements  but  rather  to  numerical  approximations  of 
measured  data.  Figure  49  and  Figure  50  display  the  semi-empirical  results  for  complex 
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permittivity  as  a  function  of  frequency  for  the  same  physical  conditions  and  soil  compositions  as 
considered  in  the  calculation  of  the  empirical  results  shown  in  Figure  42  and  Figure  43. 


COMPARISON  OF  MEASURED(SOUD  CIRCLES)  AND  CALCULATED  VALUES  OF  PERMimVITY 


Figure  45 
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Figure  49 
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Figure  50 


4.4.6.  Applied  Linear  Regression  Analysis  in  the  Frequency  Range  1-50 
MHz 

In  order  to  construct  a  comprehensive  numerical  algorithm  capable  of  evaluating  complex 
permittivity  for  arbitrary  frequency,  volumetric  water  content,  and  as  a  function  of  the  weight 
fraction  of  the  primary  solid  components  silt  sand  and  clay,  it  is  required  to  collect  sufficient  data 
which  adequately  spans  the  full  range  of  the  mixtures  understood  to  comprise  soil  solids.  The 
numerical  curve  fitting  (of  permittivity  as  a  function  of  moisture  content  and  soil  solids 
composition)  of  Dobson  [14],  provides  an  initial  model  of  the  desired  functional  dependence  but 
restricted  to  the  frequency  range  1400-18,000  MHz.  As  becomes  apparent  from  the  phase 
diagram  in  Figure  51  where  the  soil  compositions  considered  by  Dobson  are  indicated  by  solid 
black  circles  numbered  4,5,6,7  and  8  (corresponding  to  his  designations  Fields  I,  II,  II,  IV  and  V 
respectively),  the  soils  considered  are  relatively  localized  to  the  (central)  loam  and  silt  triangles 
and  cannot  be  expected  to  adequately  characterize  soil  sample  significantly  displaced  into  the 
sand  and  clay  triangles.  The  data  of  Campbell  [1 1]-[12],  in  its  consideration  of  the  primary  soil 
varieties  (Wilder)  silt,  (Fort  Edwards)  clay  and  (Hart)  sand  extends  the  data  set  to  the  primary 
extremes  of  the  phase  diagram  introducing  the  red  circles  at  each  of  the  three  vertices.  The  data 
of  Campbell  is  also  of  interest  in  its  consideration  of  frequencies  below  the  lower  limit  of  the 
Dobson  range  of  1400  MHz  to  include  measurements  at  (1,2,5)  MHz  and  (10,20,  50)  MHz. 
Experimental  data  for  even  lower  frequencies  is  found  in  the  studies  of  Sternberg  [17]  who 
explores  the  full  range  of  frequencies  0.01  -  100.00  MHz  for  Ca-Montmorollinite  (clay,)  natural 
sand  (Avra  Valley)  ,  and  for  pristine  (Brookhaven)  sand  ,  assumed  to  correspond  to  the  lower  left 
vertex  of  the  sand  triangle. 
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Figure  51 

The  data  of  Campbell  [1 1]-[12],  provided  in  graphs  as  a  function  of  volumetric  moisture  at  six 
frequencies  was  digitized  and  processed  by  linear  regression  analysis  to  extract  the  best 
coefficients  for  a  polynomial  of  third  order  in  the  volumetric  water  content.  The  results  for  the 
modeling  of  the  data  of  Campbell  (including  digitized  measurements  and  modeling  as  cubic 
polynomials)  at  each  of  the  six  frequencies  considered  by  Campbell  are  shown  in  Figure  52, 
Figure  53  and  Figure  54  for  Wilder  Silt,  Hart  Sand  and  Fort  Edward’s  clay,  respectively.  Plots  of 
the  logarithm  of  permittivity  for  otherwise  identical  physical  conditions  are  shown  in  Figure  55, 
Figure  56  and  Figure  57.  The  modeling  expresses  the  permittivity  as  polynomial  cubic  in  the 
(fractional)  moisture  content  according  to  the  expressions 

3 

£\{v,  w„.  =  c'„  ” 

n=0 


K'«v)  =  Zc 

n=0 


m.. 


(5) 
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This  formalism  is  extended  to  arbitrary  composition  by  assuming  that  the  coefficients  in  the 
above  equation  may  be  expressed  as  a  linear  combination  of  their  compositional  fi'actions  in  silt 
sand  and  clay 

|c>„=W|a>„=>|a>„=W"' |c>„  (6) 

where  the  column  vector  consists  of  a  coefficients  (|c>„)  of  the  nth  power  from  each  of  the 
three  independently  modeled  soils.  The  unknown  combinatorial  coefficients  \  a>  are  obtained 
through  inversion  of  the  composition  matrix  W.  The  composite  modeling  of  the  primary  soils  in 
the  intermediate  range  as  described  above  extends  the  lower  frequency  range  of  permittivity  as  a 
function  of  frequency  and  moisture  content  at  arbitrary  soils  solids  compositions  from  1400  MHz 
to  1  MHz. 
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Figure  57 


4.4.7.  Projected  Scaling  to  Measured  Data  in  the  Frequency  Range  .01-100 
MHZ 

Since  the  lower  frequency  data  of  Sternberg  is  only  available  for  each  soil  considered  at  a  single 
value  of  volumetric  moisture,  the  permittivity  below  1  MHz  was  modeled  as  a  scaled  projection 
of  the  behavior  at  1  MHz.  The  scaling  factor  is  determined  by  the  ratio  of  the  permittivity  at  the 
lower  frequency  and  single  specified  moisture  to  that  at  1  MHz  of  the  same  volumetric  moisture. 
The  results  of  the  scaled  projection  of  Hart  sand  and  Fort  Edward  Clay  to  Sternberg’s  data  below 
1  MHz  for  Brookhaven  Sand  and  Ca-Montmorollinite  clay  respectively  are  shown  in  Figure  58 
and  Figure  59.  Values  of  the  complex  permittivity  of  soils  of  arbitrary  texture  at  the  lower 
frequencies  are  obtained  by  first  computing  their  values  at  1  MHz  relative  to  that  of  the  known 
bounding  extremes  (sand  and  clay)  at  1  MHz.  The  absolute  values  at  the  lower  frequencies  are 
then  obtained  from  the  known  bounding  values  at  the  lower  frequency  (as  obtained  from  the 
experimental  data)  such  that  the  ratios  of  the  permittivity  of  the  arbitrary  texture  to  that  of  the 
bounding  values  (at  1  MHz  and  at  the  lower  frequency)  are  preserved. 
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Figure  59 


4.4.8.  Performance  of  the  Model  over  Full  Frequency 
4.4.8.I.  Range  for  Representative  Soil  Compositions 

Plots  of  both  real  and  imaginary  parts  of  the  eomplex  permittivity  for  each  of  the  five  soils 
considered  by  Hallikainen  and  Dobson  and  for  the  primary  compositions  considered  by 
Campbell  (using  the  semi-empirical  model  above  300  MHz)  are  shown  in  Figure  60,  Figure  61, 
Figure  62,  Figure  63,  Figure  64,  Figure  65,  Figure  66  and  Figure  67.  Where  explicit  data  is  not 
available  as  a  function  of  frequency  or  between  major  frequency  subdivisions,  the  (log)  of  the 
permittivity  is  interpolated  as  a  function  of  the  log  of  frequency  using  the  bounding  values  of 
data  of  the  relevant  subinterval.  The  data  of  Sternberg  for  Brookhaven  sand  and  for  Ca- 
Montmorollinite  is  also  shown  on  the  two  plots  corresponding  to  these  respective  varieties  of 
soils  and  are  accompanied  by  the  corresponding  calculated  values  computed  at  the  same  value  of 
volumetric  moisture  as  the  data. 
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